Exotic S = 1 spin liquid state with fermionic excitations on triangular lattice 
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Motivated by recent experiments on the material Ba3NiSb209 we consider a spin-one quantum 
antiferromagnet on a triangular lattice with the Heisenberg bilinear and biquadratic exchange inter- 
actions and a single-ion anisotropy. Using a fermionic "triplon" representation for spins, we study 
the phase diagram within mean field theory. In addition to a fully gapped spin-liquid ground state, 
we find a state where one gapless triplon mode with Fermi surface coexists with d -\- id topological 
pairing of the other triplons. Despite the existence of a Fermi surface, this ground state has fully 
gapped bulk spin excitations. Such a state has linear in temperature specific heat and constant in 
plane spin susceptibility, with an unusually high Wilson ratio. 

PACS numbers: 71.27.+a, 75.10.Jm, 75.10.Kt, 75.30.Kz 



Spin liquid (SL) is a long sought exotic state of matter 
proposed by Anderson [1], where long range magnetic or- 
der is destroyed by quantum fluctuations at zero temper- 
ature. A number of materials have been discovered which 
are promising candidates for two-dimensional 5* = 1/2 
SL state [H. More recently, possible SL materials with 

= 1 have been discussed. One example is the insulating 
spin-1 quantum magnet on a triangular lattice, NiGa2S4, 
reported by Nakatsuji et al [3]. This material motivated 
a number of theoretical papers proposing different mi- 
croscopic realizations of 5* = 1 SL (J42|- Recently high 
pressure synthesis of the two-dimensional triangular mag- 
net Ba3NiSb2 09 has produced two new phases which 
possibly realize two and three-dimensional 5* = 1 SL. In 
particular the 6II-B phase, described as a triangular lat- 
tice of Ni^+ ions, shows no magnetic ordering down to 
T = 350 mK along with linear in temperature specific 
heat (with unusually high coefficient) and constant spin 
susceptibility. The metal-like behavior of specific heat 
and spin susceptibility observed in the insulating 6H-A 
phase suggest the presence of quasiparticle excitations 
with a Fermi surface. 

Motivated by this newly discovered material, in the 
present Letter we propose a new candidate SL ground 
state with exotic physical properties. Our model system 
consists of quantum = 1 spins forming a triangular 
lattice. For simplicity, we consider only nearest neighbor 
interactions. The general form of Hamiltonian can be 
written as 



H 



(ij) 



(1) 



where we included Heisenberg exchange interaction with 
coupling J > and biquadratic exchange with cou- 
pling K. In addition we allow easy-plane or easy- axis 
type of anisotropy controlled by the parameter but 
we neglect this anisotropy in the couplings J and K 
since it is presumably small for transition metals. The 
Hamiltonian ([1]) has been considered in the literature in 
limits when the anisotropy is either zero or dominates 
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FIG. 1. Schematic representation of the ground state in dif- 
ferent limits of the Hamiltonian ([T]). White arrows represent 
average spin; arrows with discs indicate the director of the 
nematic order parameter. Details are discussed in the text. 



over other couplings, or there are longer range compet- 
ing exchange couplings. Fig. [1] summarizes known results 
for the ground state (GS) phase diagram in a schematic 
way. There are three different phases on the line of zero 
anisotropy D = in the range K = -0.4 J. . . J 

GS is 120°-degree antiferromagnet (AFM). For larger 
negative K system favors collinear ferro-nematic (FN) 
order, i.e. nematic order that does not break lattice 
translational symmetry. In this state the average spin 
vanishes (S) = 0, but full spin rotation symmetry is 
broken down to rotations around an axis specified by 
the director vector d (see Refs. lo|, El and discussion 
below). For positive K > J the ground state is de- 
scribed by aniferro-nematic (AFN) order. In this state 
director vectors di on three different sublattices are or- 
thogonal to each other (see Fig. [T]), thus breaking lattice 
translation symmetry. In the extreme case of easy-plane 
anisotropy {D ^ J, |i^|), the GS is a trivial product 
of states of = 0) on all sites, corresponding to the 
trivial single-site FN order. For large but negative 
implying extreme easy axis anisotropy, only two states 



2 



with S*^ = ±1 on each site survive. This system can be 
described by a spin- 1/2 XXZ model with ah exchange 
couphngs being antiferromagnetic if 2 J > > or with 
being frustrating and ferromagnetic \i K < {). In 
both cases there is spin density wave ordering of the z- 
component of the spin in the GS, supplemented by planar 
AFN order in former and collinear nematic order in the 
latter case [l3| . 

Physically for Ba3NiSb2 09 we may expect the ex- 
change coupling J to be the largest with J > |D|. 
Both signs of D seems plausible. Likewise it is not known 
what sign of the biquadratic exchange K is realized, even 
though negative K can be obtained from large U expan- 
sion of certain multi-orbital Hubbard model or from cou- 
pling to phonons. Therefore in what follows we study 
the phase diagram of Hamiltonian ([1]) for both signs of 
D and K but will assume |i^| < J. Except of very 
small \D\^ this is outside of the regions of known GS's 
shown in Fig. [TJ In order to to get access to the (RVB- 
like) state with fermion excitations, we use the fermion 
representation of the spin Q . After this we study result- 
ing phase diagram in the mean field approximation. 

Fermion representation. The spin operator is conve- 
niently represented via a set of three operators called 
triplons which are labeled by index a = x^y^z. In earlier 
papers 0, [sj these operators were bosons, but here we 
use fermions ^ written as a vector fi = {fix^ fiy^ fiz)^^ 



-i/7 X fi, fl-h = 1. 



(2) 



In terms of eigenstates, we used the following basis to 
represent the states of 5* = 1, \x) = i(|l) — | — l))/y/2^ 
\y) = (|1) + I - l))/v/2, 1^) = -i|0), since it facilitates 
the handling of the biquadratic term in the Hamiltonian. 
Eq. ([2]) also imposes a constraint of single occupation in 
order to exclude unphysical states from the Hilbert space. 
In the mean field theory this constraint will be relaxed 
to hold only on average. There are two possible choices 
of constraint for spin-one system: the particle represen- 
tation that we used above and the hole representation 
fj ■ fi = 2. In contrast to the case of S = 1/2, these 
are not equivalent. Nevertheless, they can be mapped 
into each other by particle-hole transformation plus a 
change of the sign of hopping. Therefore we consider 
only particle representation but do not restrict hopping 
to be positive to include the hole representation [14j. 

The chosen spin representation has U{1) redundancy 
remaining [gI, U^: one can multiply ^ by a phase fac- 
tor leaving the spin intact. In addition, in the absence 
of D there is a spin rotation symmetry, realized by the 
simultaneous rotation of the vectors fi and {f})^ - Non- 
zero anisotropy D breaks full spin rotation symmetry to 
rotation symmetry in x?/-plane supplemented by the re- 
fiection of spin along z-axis. 

The bilinear term is expressed via fermions as Si ■ Sj = 
(/7 • fj)ifi ■ fj) + Piifi ■ Pj)fj- Using the constraint /7 • 



/i = 1, the biquadratic term also can be expressed as a 
product of four fermion operators fo'], {Si • SjY = 1 — 
{fi ' fj)ifi ' fj)' Adding a Lagrange multiplier to enforce 
the single occupancy constraint (|2]) on average, we have 



H ■ 



^[Jflifi • f~])fj + - K){fJ • f])ifi ' fj) + K] 

(ij) 

+ E[Mi-/7-/^) + ^(i-/l/-)]. (3) 



Mean field results. Having expressed the Hamiltonian 
via fermion operators we study the mean field phase dia- 
gram of our model. To unambiguously decouple quartic 
fermion terms, we use the Feynman variational princi- 
ple m, [TgI which is equivalent to the trial wave functions 
approach. We define an action based on the Hamilto- 
nian daj), 6' = Jq drl^- fl^{dr - IJi)fiQ, + i^] , as weh as 
the trial quadratic action, 5, with H replaced by 

H = Y^^flTiJi + /7^i,-/] + H.C.] + ^ fluh. (4) 

The mean field parameters Tjj, Aij^ and ti are de- 
termined from the stationary points of the functional 
^[S] = {S-S)s-\ogZ, 



T^P = -J5^p{fUi.) + (J - K){fUip), 



na(3 

^ij 

af3 



-J{fi(3fja) -^{J - K) Sa(3{fiKfjK), 



^ij 

tf = T.iJ{fhfj^) - - ^)(flfj^)] - DSa^^a.^ 
(ij) 

For T = 0, we get the estimate for the ground state 
energy, E^^s. < ^g.s. = {H)s, where 



i 

We search for self-consistent solutions to the mean 
field equations that do not break any additional sym- 
metries other than T-reversal. When the full spin ro- 
tation symmetry is present, the only possible pairing 
order parameter is Aq ~ {fi ■ fj). Such pairing pre- 
serves full rotational symmetry in spin space, the result- 
ing state being a spin singlet. We call this pairing in 
odd channel, since it is possible only with odd orbital 
momentum, i.e. p, /-wave pairing. Since in Hamilto- 
nian ([T]), only in-plane rotational symmetry is present 
for D 7^ 0, the pairing in even channel with order pa- 
rameter Ae {{fi X fj)z) = {fixfjy - fiyfjx) is alloWCd. 

However, the presence of two order parameters simulta- 
neously violates the symmetry with respect to rotations 
of TT around the x or y axis. 

Both aforementioned types of pairing were considered 
by Liu et.al. [6] in a similar system, however without 
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anisotropy but with competing third nearest neighbor J. 
Their treatment of biquadratic exchange also differs from 
ours. The result of [6| was that pairing in odd channel 
always wins. Below, after establishing the mean field 
equations for each type of pairing, we identify the region 
in phase space where even-channel pairing has lower en- 
ergy than odd-channel pairing. 

Pairing in odd channel We introduce the mean field 
parameters x", and , a = x^y^z defined as 



X — {fiQ,fi-\-eia)i ^ — {fiafia)i ~ {fiafi-\-eia)' 

(7) 

The vectors ei = (1,0), €2 = (1/2,^3/2), and 63 = 
62 — ei specify link orientation. The hopping is the same 
on all links, whereas the pairings for the remaining two 
orientations are {fiafi+e2a) = A^e"?^, (fiafi+esa) = 
A^e^", where the pair angular momentum / = 1,2,3 
for p + ip, d -\- id, and /-wave pairing respectively. Spin 
rotation symmetry in the xy-plane requires = X^? 

_ ^ /\^x _ rpj^^ Hamiltonian in momen- 

tum space (modulus non-essential constant terms) can 
be rewritten as 



H 



k.a 



Xkfkafka + ^kfkaf-ka + f-kafka, (8) 



with mean field parameters 

= 27(fe) [(J - /f)x" - J(X" + X' + Xl] (9) 

+QKn°' -II- 5a,zD, 
= V(fc)[(J - K){^1 + + A^) - JA?]. (10) 

The function ^{k) is a sum over nearest neighbors, 
7(/c) = cos/c-ei+cos/c-e2+cos/c-e3. On the other hand, 
il){k) depends on the type of pairing under consideration. 
Note that p-wave pairing breaks lattice rotational sym- 
metry. Therefore we consider p + ip-wave and /-wave 
pairings: ip^ (k) = i(sin/c • ei — sink • 62 + sink • 63), 
^pip(^k) = i(sin k- ei^ e'""^^ sin k ■ €2 ^ e^^""/^ sin k - es). 
Eq. (jHj) is solved with Bogoluybov transformation act- 
ing separately on each fermion species. This results in 
the spectrum 
equations: 



and mean field 



N 



1 



2Ei 



k 



=-y- 

N ^0 



2E^ 



2E^ 



(11a) 



(lib) 
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supplemented by the constraint equation {fj - fi) = I. 

Pairing in even channel Hoppings are defined as 
in (jlj), whereas pairing is A^^ = l/2{fixfi+eiy - 
fiyfi+eix)' The Hamiltonian is: 

^ ~ ^^^Xkfkaff^(^ + ^k fkxf-ky + ^k f-kyfkxi 
k,a 



with Xk given by Eq. (j9j), and A^^ = 2Jilj{k)A^y. Note, 
that the fz band is unpaired and retains its Fermi sur- 
face. We consider 5-wave and d + id-wave pairings {d- 
wave violates lattice symmetry and higher orbital mo- 
mentum pairing requires inclusion of further neighbors). 
For the case of s-wave pairing, the function ilj^{k) = j{k). 
For d + id- wave pairing we have ip^^^{k) = cosk • ei + 
cos k-e2-\- e~^^^/^ cos k es. The Bogolyubov spec- 



p2i7r/3 



trum is El 



xi 



. Self- 



consistent mean field equations for x and ^/-components 
are given by Eq. (pTj) with the new expressions for the 
spectrum and gap functions. For the z component we 
have 

X' = J^Y1 ^7(^)^F(Xfe). = ^ ^f{xI)' 

k 



Our mean field approach includes on-site FN order au- 
tomatically. The on-site nematic order is described by 
the order parameter tensor, = 1/2{S''S^ ^ S"") - 
2/36^^. For a single site with /S = 1 all states with 
zero average spin (S) = can be characterized by the 
unit director vector d [lo(, in the basis defined earlier, 
\d) = dx\x) -\- dy\y) -\- dz\z). For this state Q"^ is ex- 
pressed via d as Q^^ = l/3dai3 — dcxdf3. For example, 
d\\z corresponds to the state \S^ = 0), and the nematic 
order is diagonal, Q^^ = diag(l/3, 1/3, — 2/3). In our 
model we also have states with vanishing spin order and 
diagonal on-site nematic order. However, since our GS is 
RVB-like with long-range entanglement, Q^^ cannot be 
described by the above simple form. We have to intro- 
duce the magnitude Q^^ = q{l/3Sai3 — d^d^). Calcu- 
lating the nematic order parameter tensor in our model 
we have = Jq,/3[1/3 — n"], where is the average 
occupation of corresponding fermion. Since = n^, we 
have nematic order with d||z, with a magnitude given by 
q = — n^^ varying from 1 for = 1 (state = 0)) to 
— 1/2 for = 0. Non-zero anisotropy D causes 
to be different from 1/3, and therefore directly couples 
to FN order along the z-axis. 

Having studied the energies of all aforementioned 
states using Eq. (|6j), we found that the main competition 
is between states with p -\- ip and d + id- wave pairings, 
all other states being higher in energy. As one increases 
the effective coupling for the odd-channel pairing de- 
creases, whereas for even pairing it remains the same. Fi- 
nally for K ^ 0.45J singlet pairing wins. The resulting 
phase diagram is shown in Fig. [21 The boundary between 
two states appears to be weakly dependent on D. 

Physical properties of d -\- id state. The d + id state 
breaks time-reversal symmetry. The chiral order param- 
eter associated with this broken symmetry {Si • {Si-^^i x 
Si+es)) oc x^|Ag^p is proportional to the magnitude of 
the pairing gap squared. In addition, pairing with d -h id 
gap symmetry in two dimensions is topological [l^, re- 
sulting in the existence of a pair of zero-energy edge 
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FIG. 2. The phase boundary between SL GS's with p-\-\p and 
d -\-\d pairing, (a) The spin susceptibihty Xxx in the d -\- \d 
phase as a function oi D / J for K/ J — 0.55. The susceptibihty 
is normahzed by the average density of states, v — (vx^Vz) l^-, 
where is calculated without the gap. (b) Gapped (dashed 
red line) and ungapped (blue line) Fermi surfaces of x, ^, and 
^-fermions for KjJ = 0.55, DjJ = 0.8. 



scale. These results tell us that the Fermi surface as- 
sociated with [see Fig. [2] (b)] should be viewed very 
differently than the spinon Fermi surface in the 5* = 1/2 
SL which carries spin- 1/2 quantum numbers and leads 
to gapless spin-1 excitations. In our case 5*^ = 1 exci- 
tations are gapped even though the static spin suscep- 
tibility Xxx^Xyy 7^ and the specific heat has linear T 
dependence. 

Finally, we discuss experiments that could confirm the 
proposed ground state. Measurement of the spin suscep- 
tibility for single crystal or oriented powder samples is 
of great interest in order to test our prediction of strong 
anisotropy. We also predict an exponentially activated 
behavior for l/(TiT) which may be surprising in view of 
the linear T behavior of the specific heat. 

We thank Luis Balicas for bringing Ref. [H to our at- 
tention. We acknowledge useful discussions with Samuel 
Bieri. T.S. is supported by grant NSF-DMR 6922955. 
P.A.L. is supported by NSF-DMR 1104498. 



modes at the boundaries. The physics of these modes 
will be discussed elsewhere. 

The combination of gapless excitations with topolog- 
ical pairing gives rise to a number of unusual physical 
properties, that may explain the results of the recent ex- 
periment [8]. Due to ungapped fz excitations the spe- 
cific heat depends linearly on temperature near T = 0, 
C = 7r'^k%i'zT/3, where Vz is the density of states of 
fiz at the Fermi surface. Due to Higgs mechanism the 
gauge field is massive and does not modify the linear in 
T behavior of the specific heat. The spin susceptibility 
exhibits more exotic behavior: due to the pairing of x 
and ?/-fermions the zz-component Xzz = 0. On the other 
hand, Xxx is finite and depends on the anisotropy D. 
For D smaller than the gap, Xxx = XxxKi^bq)'^ ~ ^z, 
and approaches a factor two larger value Xxx ^ 2^^^, 
when D is much larger than the gap. This difference 
by factor 2 is approximate, valid in the limit of con- 
stant gap and density of states. The behavior of Xxx 
is shown in Fig. [2] (a). We calculate Wilson ratio de- 
fined as Rw = (47r^/c^)/(3^^/i^)(x^)/C, and obtain 
Rw = 8/3 ~ 2.66 for the case of small anisotropy, and 
Rw 16/3 ^ 5.33 for large anisotropy. Note that we 
take the average susceptibility X = 2/ 3xxx to account for 
the polycrystalhne nature of the sample. The latter value 
gives surprisingly good agreement with the Wilson ratio 
observed experimentally, Rw ^5.63. We also calculated 
the imaginary part of the spin susceptibility. Since two 
out of three fermions are gapped, lmxaa{^^<l) vanishes 
for temperatures and frequencies smaller than the gap 
for all a. This implies the NMR relaxation l/(TiT) is 
exponentially small for temperatures below the pairing 
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